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SIMULATIONS  OF  A  MISSILL  WITH  LIQUID 
PAYLOAD  UNDER  LAUNCH  CONDITIONS 


INTRODUCTION 

This  project  was  concerned  with  the  numerical  simulation  of  the  flow  of 
a  liquid  in  the  annular  region  between  two  coaxial  circular  cylinders. 
The  Investigations  related  to  the  behavior  of  a  missile  with  a  liquid 
payload  during  and  immediately  following  launch. 

The  project  comprised  two  Parts.  In  Part  1  solutions  were  computed  of 
the  transient  Navier- Stokes  equations  for  a  number  of  different  liquids, 
both  Newtonian  and  non -Newtonian,  with  the  object  of  determining  whether 
or  not  the  liquid  became  fully  spun  up  over  a  given  Interval  of  time.  In 
Part  II  the  steady- state  Navler-Stokes  equations  relative  to  an  aerobai- 
llstlc  frame  were  solved  with  the  aim  of  determining  the  despin  moment 
on  a  spinning  md  nutating  missile  in  a  fully-spun-up  mode. 

Detailed  descriptions  of  the  two  problems  will  be  given  below.  In  both 
Parts  the  simulations  were  performed  with  the  finite  element  software 
package  FIDAP,  Most  of  the  computations  were  done  on  an  HP  9000  mini¬ 
computer;  some  were  done  on  a  Mlcrovax  workstation. 
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Part  1 


SPIN  UP 


Computations  wero  performed  to  determine  the  spin  up  chnrnctcrist  I cs  t)f 
a  number  of  different  liquids,  both  Newtonian  and  non-Ncwtonl an ,  undor 
conditions  simulating  the  launch  of  a  missile  with  a  liquid  payload. 

The  essential  features  of  the  process  studied  were; 


the  time  interval  over  which  spin  up  took  place  was  from  0.01  see 
to  0.12  sec; 

the  spin  rate  was  accelerated  over  this  time  intervMl  in  a  pre¬ 
scribed  manner  (detailed  below)  from  zero  to  11.8  rov/soc. 


the  axial  acceleration  varied  over  the  same  time  interval  in  a 


prescribed  manner  (detailed  below)  ,  with  a  value  of  b81)  ir/sec^  at 


the  end  of  the  interval. 


The  liquid  was  taken  to  occupy  the  region  between  two  coaxial  circular 
cylinders  of  radii 


iJj  -  2.625"  .  •  4.5" 


respectively.  The  boundaries  of  the  annular  region  were  rigid  and  the 
length  of  the  cylinders  was  66".  However,  the  liquid  filled  only  70«  of 
the  region,  so  that  the  effective  (mean)  length  of  the  liquid  column  w.u; 


d  -  46. 2 ■ 


The  equations  governing  the  motion  of  the  liqviid  are  the  Navier- Stoke.s 
equations  and  the  continuity  equation,  which  can  be  written  a.s 


/’(y\*  -t  y*-vV)  -  -'7V  +  v*.7*  t  pf: 


V  V  -  0 


a 


.'k  .-''I.-* "Jir^  V  V  V  v  '.»'.’J<  V- 


vhcie  />  is  the  density,  y*  Is  velocity,  p*  is  pressure,  is  extra 

stress  and  £*  Is  a  body  force  vector  which  depends  on  the  time  f * .  As£:o- 
clated  with  those  equations  Is  a  constitutive  relation  between  the 
stress  tensor  and  the  str&ln-rate  tc-soi ; 

(3) 


(^) 

For  Newtonian  liquids  (3)  Is  Just  a  linear  telarlon;  for  the  non-Ncwto 
nlan  liquids  studied  was  given  by  empirical  data. 


I*  "  ^  sC  (D 

where  Pp  Is  zero- shear -rate  viscosity  and  where 
jf  -  h(V  j-  +  V  i  ^) . 


Nondlniensloriallzation 

It  is  convenient  to  express  the  governing  equations  and  associated  i;uaii- 
titles  In  nondlni^nsional  form.  All  lengths  are  scaled  with  re.spcct  to 
the  outer  radius  R^,  time  with  respect  to  the  final  angular  velocity 
w  (-11.8  rev/sec),  and  velocitle.s  wit),  respect  to  the  velocity  R^w.  Thus 
we  write 

ii*  -  3?jU)V,  r*  -  w'*t,  c  -  R^c.  p'  -  pC^juO'^p, 

t  -  ^<^1.  .  2  '  “'2 

Then  the  governing  equations  (1)  •  (4)  become 

■»  Ji-Vi:  -  -Vp  +  ]^  V  T  +  t(f)  O) 

V-v  -  0  (6) 


with 


I  -  £(2) 


(/) 


and 


2  —  h(Vv  ^  Vv^) 


(8) 


9 


The  paramctr-r  appearing  In  equation  (b)  Is  the  Reynolds  number  based  on 
zero-shear-rate  viscosity. 


Re  - 


(9) 


We  now  specify  In  more  detail  certain  scalings  used  in  the  calculations. 
The  final  angular  velocity  is 

o  -  11.8  revs/sec  -  74.14  rads/sec 


so  that 


ui’^  -  .013b  sec. 


The  time  of  spin  up  is  from  f*-  O.Ol  sec  to  t*-  0.12  sec,  and  the  spin 
rate  increases  linearly  from  (*-  0.01  sec  to  t*—  0.08  sec,  and  then  is 
constant.  We  therefore  define  a  dimensionless  time  t  by  the  formula 

f‘  -  .01  +  ,0135c 


or,  equi'  .UfMtly, 

The  relation  between  dimensional  and  dimensionless  spin  up  times  is 
therefore  a.s  shown  in  the  following  table 


c" 

f 

.01 

0.0000 

.08 

5 .1900 

.  12 

8.1554 

Thus  the  total  dinicns  ionless  spin  up  time  is  r  -  8.1354. 


'■j'  1 

I 

.-vl 

VO 


^,V- 

.  e.  > 
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Wf-  can  also  conveniencly  delino  n  tliimii.-j  ion  less  angular  velocity  tJ  which 
changes  In  time,  as  Indicated  abov'e,  as  follows: 


n  -  n(;) 


t/5.19  for  o  <  r  <  5.19 

1  for  5.19  <  t  <  8.1554 


(10) 


It  is  convenient  also  to  non-dimonsionalize  the  time -dependent  body 
force.  This  increases  until  t*“  0.08  and  then  experiences  a  small 
decrease.  At  the  final  instant  t.*’-  0  12  it  is  given  that 

F*  ,  -  585.50  m/sec^ 
irn«l  ' 

When  this  Is  scaled  with  respect  to  we  obtain 

f'r,  ^  -  K  ,  “  0.9319 

zinal  iinal  Z 

Hence  we  write  the  dimensionless  acceleration  over  the  time  interval  of 
Interest  as 


f(c)  -  -  0.9319  /(r) 


where 


/<c) 


.1115  +  .17901 
1.1108  -  0.0136r 


0  <:  c  <  5,19 
5.19  <  c  <  8.1554 


(11) 


(12) 


The  values  at  different  times  rt-proseiiLed  by  (11)  and  (12)  are  indicated 
in  the  following  table; 


^  1 

f(0 

* 

) 

m 

i 

F*(t*) 

0 

.111') 

-  .103'^ 

.01 

-  45.28 

5.19 

1 .0403 

-  96-)  3 

.08 

-  609.1 

8.1554 

1  .  OUOU 

-  .93)9 

.  12 

-  585.5 

The  s  imulal  iorii,  were  perforinco  t<.i  .i  ivi.il  of  1:)  ditlcienl  liquids,  hive 
were  Newtonian  and  ten  weif  non  Nv  w  ton  i.iu  hor  i  dent  i  i  i  ca  t  i  on  purpost  .s 


the  ten  non-Newtonian  liquids  were  labeled  i iquid  1,  liquid  2 . liquid 

10.  The  Newtonian  liquids  were  labeled  liquid  IN,  liquid  2N ,  liquid  5W, 
liquid  6S  and  liquid  SN .  This  was  because  of  equality  between  the  vis¬ 
cosity  of  liquid  IN  with  the  zero-shear-rate  viscosity  of  liquid  1,  etc. 
These  viscosities  and  the  corresponding  Reynolds  numbers  are  shown  in 
Table  1.  The  dependence  of  the  non-Newtonian  liquids'  viscosities  on 
shear  rate  was  determined  from  empirical  data  supplied. 


liquid 

n 

Mo 

poise 

Re 

1 

1000 

9.688 

IN 

1000 

9.688 

2 

600 

16.15 

2N 

600 

16.15 

3 

380 

25.50 

A 

200 

68.46 

5 

100 

96.88 

5N 

100 

96.88 

6 

60 

161.5 

6N 

60 

161.5 

7 

20 

686.6 

8 

10 

968.8 

8N 

10 

968.8 

9 

6 

1615 

10 

2 

6866 

TABLE  1.  Liquids  and  zcro-shcar-ratc  viscosities 
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The  Initial  conditions  for  the  simulation  were  that  the  liquid  was  at 
rest  at  (dinens ionless)  time  t  -  0. 

The  boundary  conditions  are  applied  to  all  the  surfaces.  The  inner  and 
outer  curved  boundaries  are  at  dimensionless  radial  positions  0.5833  and 
1.0  respectively.  The  lower,  solid  end  boundary  of  the  cylinder  is  at 
z  ■*  0  in  a  cylindrical  polar  coordinate  system,  while  the  mean  position 
of  the  upper  free  surface  is  at  z  10.27. 

The  boundary  conditions  are: 

1.  On  the  inner  curved  surface  r  “  0.5833,  0  <  z  <  10.27, 
u  -  w  -  0  ,  V  -  0.5833  Q(C) 

where  0(C)  is  given  by  equation  (10);  (u,v,w)  are  the  velocity 

components  in  the  radial,  azimuthal  and  axial  directions  respec¬ 
tively. 

il.  On  the  outer  curved  surface  r*=1.0.  0<z<  10.27, 
u  -  w  -  0  ,  V  -  l.Q(c) 
ill.  On  the  base  z  »  0, 

u  “  w  -  0  ,  V  -  rfi(£)  ,  0.5833  <  r  <  1.0. 

IV.  On  the  upper  free  surface  it  is  assumed  that  the  surface  tension 
is  zero.  The  boundary  conditions  on  this  surface  are  then  that  the 
normal  velocity  is  zero  and  that  the  tangential  stress  is  zero. 
These  translate  into 


V  -  ^  =  0 

dzdz 


on  Che  mean  position  z  -  10.27 


The  problem  was  solved  as  a  transient  problem  over 
interval 


a 


WWW 


M  »■  m.  J I  k  H  L  u.  LT.V’.y.S,'  WA'-jy.V.V.V.VAV^T 


7^  j.^^r^'^yvyf^,yyvj,yjr^j,-2 


0  <  t  <  8.1554.  Fixed  time  increments  were  used  of  length  0.02039,  .so 
that  40  time  steps  were  required  to  complete  the  simulation.  Each  run 
took  approximately  8  hours  of  c.p.u.  time  on  an  HP  9000  computer,  except 
for  liquids  9  and  10,  where  the  run  time  became  excessively  long  because 
of  the  need  to  design  an  extremely  fine  mesh. 

The  results  are  presented  in  the  accompanying  graphs.  Figures  1-13  show 
the  behavior  of  the  azimuthal  component  of  velocity  at  one  particular 
point  in  the  liquid  as  a  iL.nction  of  time.  This  point  was  chosen 
(arbitrarily)  to  have  a  radial  coordinate  r  -  0.756.  The  velocity  at 
this  point  was  zero  at  t  -  0,  and  if  the  liquid  is  fully  spun  up  into 
rigid  body  motion,  it  should  have  an  azimuthal  velocity  v  -  0.756  at  the 
end  of  the  time  Interval  t  -  8.155.  The  figures  show  that  the  Newtonian 
liquids  IN  and  2N  are  indeed  fully  spun  up,  while  none  of  the  others  is 
completely  spun  up.  The  non-Newtonian  liquids  1  and  2  are  nearly  spun 
up,  to  within  a  few  percent.  The  velocity  at  this  point  which  is 
actually  reached  is  a  measure  of  the  extent  to  which  complete  spin  up 
has  been  achieved.  Using  the  ratio  of  velocity  attained  to  the  spin  up 
value  (.756)  as  a  simple  measure,  we  can  refer  to  percentages  as  shown 
in  Table  2.  Liquids  9  and  10  are  so  far  from  being  spun  up  that  the 
solutions  were  not  calculated;  they  are  in  fact  less  than  1%  spun  up. 
Table  2  shows  the  dramatic  decrease  in  the  extent  of  spin  up  as 
viscosity  decreases,  as  well  as  the  significant  effects  of  shear 
thinning. 
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Part  II 

DESPIN  MOMENT 


Formulation 

Computations  were  perfomed  to  determine  the  despLn  moment  on  the  con¬ 
tainer  described  In  Parc  I  on  the  assumption  that  the  liquid  had  become 
fully  spun  up.  and  that  the  system  was  moving  under  spin  and  nutation  at 
a  prescribed  coning  angle.  The  angular  velocity  of  the  spin  Is  denoted  - 
and  the  angular  velocity  of  Che  nutation  is  denoted  11.  The  axis  of  the 
cylinder  Is  inclined  to  the  vertical  at  the  coning  angle  8.  Ttie  two  axes 
of  rotation  Interaect  at  the  midpoint  of  the  base  of  the  cylinder. 

The  problem  was  solved  by  formulating  the  governing  equations  in  a 
nutating  frame  of  reference,  the  so-called  aeroballlstlc  frame.  In  this 
coordinate  system,  the  2>axls  coincides  with  the  axis  of  the  cylinder, 
the  x-axls  lies  in  the  plane  containing  the  angular  velocity  vectors  u 
and  Q,  and  the  y-axls  is  perpendicular  to  this  frame.  This  constitutes  a 
right-handed  cartesian  system. 

The  governing  equations  were  written  in  dimensionless  form,  with 
lengths,  velocity,  pressure,  stress,  strain  rate  and  body  force  made 
non-dlroenslonal  according  to  the  same  scalings  as  were  used  in  Part  I. 
The  angular  velocity  vector  Q  of  nutation  is  non-dlmensionallzed  with 
respect  Co  its  magnitude  Q.  The  governing  equations  for  steady  motion  In 
the  aeroballlstlc  frame  are  then  found  to  be 

yVy  +  2f]Qxy  +  ij^UxCQxr)  -  -Vp  +  +  E  (13) 

7-y  -  0  (1<1) 


where  Re  is  the  Reynolds  number  defined  by  equation  (9)  and  r]  is  the 
spin  ratio 


The  dimensionless  angular  velocity  fi  with  respect  to  the  cartesian 
aeroballistic  frame  is 


Q - + 


(16) 


where 


a  —  sin  9  (17) 

is  the  third  principal  parameter  of  the  problem. 

The  boundary  conditions  imposed  for  the  solution  of  the  problem  were: 

i)  the  normal  velocity  component  is  zero  at  all  rigid  boundaries; 

ii)  the  tangential  velocity  at  rigid  boundaries  is  equal  to  the  rota¬ 
tional  velocity  of  that  boundary; 

lii)  the  normal  velocity  and  the  tangential  stress  are  zero  at  the  free 
surface  of  the  liquid. 


Computations  and  Results 

Solutions  were  computed  for  six  different  cases:  three  of  these  related 
to  Newtonian  liquid  IN  of  Part  1,  and  three  to  non-Newtonian  liquid  1. 
These  liquids  have  (zero-shear- rate)  viscosity  of  about  1000  poise  and, 
a.s  shown  in  Part  1,  were  both  fully  spun  up.  The  various  cases 
considered  reflected  the  spin  rate  and  coning  rate  at  different  instants 
during  flight. 

The  quantity  of  interest  was  the  despin  moment  denoted  in  dimension¬ 
less  form  and  rt*  in  dimensional  form.  The  relation  between  the  two,  in 
accordance  with  the  scaling  of  Part  1,  is 
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Once  the  velocity  fii'Ui  ha?.  h*on  computed  using  FIDAP,  the  despin  moment 
is  calculated  with  the  tild  ut  a  simple  subruutiite. 

The  results  of  the  computations  are  shown  in  Table  J . 

Newtonian  Liquid 


Case 

# 

rpm 

n 

rpm 

n 

spin 

ratio 

Re 

M* 

z 

ft-lbs 

1 

640 

54 

.0844 

8.75 

.000324 

.0211 

2 

255 

230 

.9020 

3.49 

.009453 

.0975 

3 

272 

136 

.5000 

3.72 

.003188 

.0374 

Non-Newtonian  Liquid 


Case 

# 

w 

rpm 

Q 

rpm 

V 

spin 

ratio 

Re 

It 

ft-lbs 

1 

640 

54 

.0844 

8.75 

.001384 

.0899 

2 

255 

230 

.9020 

3.49 

.044821 

.4623 

3 

272 

136 

.5000 

3.72 

.011733 

.1377 

TABLE  3.  Despin  Moments 

It  is  clear  that  the  non-Newtonian  liquid  exerts  a  despin  moment  about 
four  times  as  great  as  tlie  Newtonian  liquid  with  the  same  zero-shear- 
rate  viscosity,  all  other  parainei.ers  being  equal.  This  can  presumably  be 
attributed  to  the  shear  t’.  inning  effects  present  in  the  non-Newtonian 
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liquid  and  the  corresponding  change  in  Reynold*;  number.  However,  as  can 
be  seen  from  Table  3,  Reynolds  number  is  not  the  only  parameter  con¬ 
tributing  to  variation  in  despin  moment,  so  that  the  nonlinear  modifica¬ 
tion  of  the  velocity  field  in  the  non-Newtonian  case  must  be  regarded  as 
a  significant  contributing  factor  to  the  outcome. 


